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PATENT APPLICATION 



Docket No. : D440 

Inventor (s): Randal K. Douglas & Anthony S. Abbott 

Title: Spread Spectrum Receiver Kalman Filter 
Residual Estimator Method 

SPECIFICATION 

Statement of Government Interest 

The invention was made with Government support under 
contract No. F04701-00-C-0009 by the Department of the Air 
Force. The Government has certain rights in the invention. 

Reference to Related Application 

The present application is related to applicant's 
copending application entitled Global Positioning System and 
Inertial Measuring Unit Ultratight Coupling Method, S/N: 
09/396,105 filed 09/14/99, now on the military critical 
technology list, for official use only. 
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Field of the Invention 

The invention relates to the field of communications and 
navigation systems. More particularly, the present invention 
relates to communication systems coupled to inertial navigation 
systems having robust Kalman prefiltering generating residual 
estimates used for improved signal tracking. 
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Background of the Invention 



Two fundamental components of a GPS signal are a 
pseudorandom noise code and a carrier. Navigation data are 
derived from this signal by correlating it with code and 
carrier replicas generated by the GPS receiver. A code tracking 
loop attempts to drive the code correlation error, the 
difference between the receiver generated replica and the 
received signal, to zero by advancing or retarding the replica 
code rate. Receiver to satellite range, or pseudorange, is 
computed from the code phase because the code phase is 
determined by the time required for the signal to travel from 
the satellite to the receiver. When the code correlation error 
is small the receiver is said to be in code lock. 

The carrier loop tracks either the carrier frequency or 
phase by beating a replicated carrier against the received 
signal to generate an error harmonic. The measured frequency of 
the error is fed through a filter to form a correction to the 
carrier replica oscillator. The receiver is said to be in 
carrier lock when the error harmonic frequency is small, 
typically less than ten Hertz. Receiver to satellite range 
rate, or pseudorange rate, is computed from the replicated 
carrier frequency offset, or Doppler shift, from the nominal 
carrier frequency. A 50.0 Hertz navigation message superimposed 

on all GPS signals admits a 180° degree phase ambiguity in the 
carrier. A frequency locked loop or a phase locked loop 
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implemented as a Costas loop allows tracking with no loss of 
lock across data bit changes- Both code and carrier tracking 
loops adjust the phase and frequency of the generated replicas 
to account for many factors, which include receiver and 
satellite motion, receiver clock drift, and ionospheric and 
tropospheric signal delays. 

Nonzero tracking errors appear as unmodeled range and 
range rate errors that corrupt the navigation solution. With 
navigation updates computed nominally once per second, only the 
most recent data taken prior to the measurement time are 
included in measurements provided to the navigation filter. 
Most of the correlation data taken over the one-second interval 
are used to keep the tracking loop in lock. This is determined 
by the tracking loop bandwidth, which, if set too low, will not 
track the signal through receiver acceleration or turns. 

In traditional tracking loops, there is the possibility 
that the errors in the loops are time correlated because there 
is a delay in the loops as the loops drive the error signal to 
zero. Depending on the order of the tracking loop, there may 
actually be a steady state tracking error when the order of the 
dynamic motion exceeds the order of the tracking loop. This is 
commonly called the dynamic stress error. Even with inertial 
aiding to the tracking loop, the tracking error can be time 
correlated because of time correlated errors in the aiding 
data. This possibility has limited tracking performance 
because the tracking loop measurement errors can manifest 
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themselves as errors in the state vector after processing, 
which in turn can affect future aiding data and possibly cause 
the IMU navigation and GPS tracking system to generate a wrong 
but consistent navigation solution. With normal tracking 
loops, the GPS measurements used by the navigation filter are 
deliberately spaced out in time to avoid the unwanted temporal 
correlations of the tracking loops. 

In conventional tracking loops, the measurement noise 
variance used by the navigation filter is usually estimated 
from the tracking state and the carrier to noise ratio. There 
is no guarantee that these values are consistent with the 
actually realized measurement errors. Hence, a large 
integration navigation filter may process measurements with an 
erroneous uncertainty variance of the measurement noise. The 
navigation filter can be implemented as a Kalman filter. 

In conventional tracking loops , the measurement of 
pseudorange acceleration is not considered because there is no 
apparent way to obtain a low noise estimate of the pseudorange 
acceleration from the tracking loops. Conventional tracking 
loops do not have an ability to determine the pseudorange 
acceleration. The use of a pseudorange acceleration 
measurement has many potential benefits in the observability of 
user clock and IMU error instability. With a pseudorange 
acceleration measurement, there is a more direct measurement of 
higher order errors with increased observability. Normally, 
these errors are observed through dynamic coupling into 
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pseudorange and pseudorange rate measurements. Hence, without 
the pseudorange acceleration measurements, the observability of 
slock and IMU instabilities is limited. This limitation can be 
reduced using highly accurate user clocks and inertial 
measuring units. 

Advanced tracking loop designs have been proposed to 
alleviate some of the loss of lock and reacquisition problems 
by tightly coupling IMU data with the carrier and code tracking 
loops and centralizing the navigation and correlation 
functions. With the navigation filter essentially forming part 
of the tracking loop, the computational burden becomes very 
high because the navigation filter that now has many states 
must also be updated at tens of Hertz rather than the 
conventional one Hertz . 



With unlimited navigation filter processor throughput, the 
navigation filter could be operated at extremely high rates 
with no loss of optimality due to time correlated tracking loop 
errors. This high rate operation is similar to what has been 
called a vector delay lock loop. A vector delay lock loop is a 
method to track all satellites in view simultaneously with one 
common algorithm. The vector delay lock loop broadens the 
normal aided and unaided tracking loop design approach to both 
code and carrier tracking on all in view satellites. The entire 
algorithm must run at very high processing rates because there 
is no provision for federated processing. Unfortunately, 
current processor throughput cannot support Kalrnan filter rates 
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of several tens of Hertz with the large state vectors required 
for GPS inertial navigation. Large Kalman filters may be 
decomposed into one or more federated Kalman filters within a 
Kalman filter processing architecture. For example, a large 
Kalman filter can be decomposed into two partitions including a 
large integration Kalman filter and a high rate optimal 
prefilter that are more compatible with modern processing speed 
requirements. The fundamental principle is to decompose the 
complete formulation into suitable partitions such that the 
important bandwidths and models are appropriate for each 
partition. In IMU aided GPS sets designed to date, large 
navigation filters have not been decomposed into suitable 
partitions to allow optimal processing of the raw GPS samples 
with the IMU samples. 

Another approach to aiding the code and carrier 
correlation process is to design a tracking filter that 
provides estimates of carrier phase tracking error, rate and 
acceleration, satellite to user range error, and signal 
amplitude. The filtered estimates, which are processed at a 
high rate, are provided to the navigation filter where they are 
processed at a much lower rate to form the navigation solution. 
An extended Kalman filter of this design has fifth order or 
higher order dynamics and an associated matrix Riccati equation 
that must be integrated with each correlator datum. The 
computational burden is extraordinarily high. 
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Tightly coupled navigation solutions require estimation of 
the residual error in real time for computation of carrier 
phase error and code phase error during carrier demodulation 
and autocorrelation. Computational requirements of Kalman 
filter residual estimation suffer from computation complexity. 
These and other disadvantages are solved or reduced using the 
invention . 
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Summary of the Invention 

An object of the invention is to provide a method for 
improved code phase tracking in spread spectrum communication 
systems . 

Another object of the invention is to provide a method for 
improved code phase tracking and carrier phase tracking in 
spread spectrum communication systems. 

Yet another object of the invention is to provide a method 
for improved code phase tracking and carrier phase tracking in 
spread spectrum communication systems using a tightly coupled 
residual estimating Kalman filter. 

Still another object of the invention is to provide a 
method for improved code phase tracking and carrier phase 
tracking in GPS spread spectrum communications using a tightly 
coupled residual estimating Kalman filter in a GPS navigation 
system. 

A further object of the invention is to provide a method 
for generating code phase tracking errors and carrier phase 
tracking errors respectively for code phase tracking and 
carrier phase tracking of GPS spread spectrum communication 
signals using a residual estimating Kalman prefilter operating 
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on Ricatti matrices for improved autocorrelation locking in a 
GPS navigation system. 

The present invention is a communication method directed 
to carrier phase and code phase tracking using Kalman filter 
residual estimation well suited for spread spectrum 
communication systems such as GPS navigation systems. Inphase 
and quadrature correlation data are processed to provide 
estimates of carrier phase tracking error and rate and 
acceleration, as well as code phase tracking error, and signal 
amplitude. The filter state is a tracking residual applicable 
to navigation solution correction in ultratight GPS coupling 
with inertial measurements. The tracking residual estimation 
drives the code and carrier replica oscillators in tightly 
coupled correlation loops providing adjusted early and late 
code replicas and adjusted demodulation carriers for closed- 
loop code and carrier tracking. The close-loop code and carrier 
tracking can also be applied to a weak lock navigation systems 
also using the tracking residual estimation to the code and 
carrier replica oscillators . 

The method improves interference robustness and navigation 
accuracy. Improved interference robustness is realized first, 
through the use of a code tracking and carrier tracking Kalman 
filter with a Riccati matrix computation process wherein 
tracking is achieved without the need for consistent signal 
lock in the usual sense. Second, all correlator early and late 
inphase and quadrature input to the Kalman filter contribute to 
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form the measurement residual. Third, since correlator early 
and late inphase and quadrature do not need to be processed in 
real time to keep tracking loops in lock, application of more 
advanced algorithms in low signal to noise conditions becomes 
feasible. 

Navigation accuracy is improved by real time computation 
of a residual variance using the Ricatti matrix. The real time 
computation allows a navigation processor to rapidly compute 
navigation solutions in the presence of noisy signals. Range 
acceleration availability enables rapid navigation solution 
convergence with improved observability of signal 
instabilities. Using correlation data taken over a 
computational time interval provides accurate tracking error 
estimates . 

A tracking filter for the code and carrier correlation 
process of a GPS receiver provides estimates of carrier phase 
tracking error, rate and acceleration, range error, and signal 
amplitude. Both inphase and quadrature measurements are 
processed by the Kalman filter residual estimator. The filter 
dynamics and associated Riccati equations that are independent 
of the phase error enable implementation of steady state point 
designs. The Kalman filter state is a tracking residual 
applicable as a navigation solution correction well suited for 
ultra tight GPS and IMU coupling. The tracking residual is used 
to drive the code and carrier oscillators for adjusting early 
and late code replicas and demodulation carrier signals for 
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improved interference robustness and navigation accuracy. The 
filter state is an observation measurement residual that can be 
directly used by the navigation filter as a navigation solution 
correction. The error dynamics are independent of carrier phase 
tracking error. Harmonic terms in the phase error are not 
present in the filter Riccati equation nor in the error 
dynamics for enabling implementation of a steady state design. 
The real time Kalman filter processing uses filter state 
estimates to directly command both code and carrier correlation 
processes. These and other advantages will become more apparent 
from the following detailed description of the preferred 
embodiment. 
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Brief Description of the Drawings 

Figure 1 is a block diagram of a GPS inertial navigation 
system. 

Figure 2A is a block diagram of a weak lock GPS navigation 
processor. 

Figure 2B is a block diagram of an ultratight GPS 
navigation processor. 

Figure 3 is a flow diagram of a GPS signal residual 
estimator process . 
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Detailed Description of the Preferred Embodiment 

An embodiment of the invention is described with reference 
to the figures using reference designations as shown in the 
figures. Referring to Figure 1, a global positioning system 
(GPS) inertial navigation system includes a navigation system 
10 receiving GPS signals and inertial measurement unit samples 
from a sensor assembly 12 and provides position and velocity 
data to a control and display unit 14 . The navigation system 
10 functions as a GPS inertial navigation system tracking GPS 
signals from a plurality of in view satellites, not shown. The 
sensor assembly 12 includes an antenna 16 for receiving and 
providing received GPS signals and includes an inertial 
measurement unit (IMU) 18 providing the IMU sample signals, 
both of which are communicated to the navigation system 10. 
The navigation system 10 includes a downconverter 20 for 
frequency downconversion of the received GPS signals using a 

reference oscillator 24 providing an f Q internal frequency 
reference, and includes an analog to digital (A/D) converter 22 
communicating digitized GPS samples to a navigation processor 
26. The navigation processor 26 also receives the IMU samples 
and provides the position and velocity data to the control and 
display unit 14 . 

Referring to Figures 1, 2A and 2B, and more particularly 
to Figure 2A, the navigation system 10 may be a weak lock 
navigation system that includes a navigation processor 26 
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receiving the digitized GPS samples 30 that are communicated to 
mixers 32 and 34 for providing inphase (I) and quadrature (Q) 
signals respectively using cos and sin demodulation signals 
from a carrier numerically controlled oscillator (NCO) 36 using 

the reference oscillator f 0 24 . The cos and sin demodulation 
signals are replica carrier signals for carrier demodulating 
the received GPS samples 30. The I and Q sample signals are 
received by a correlator 40 receiving early (E) and late (L) 
code replica signals from a chipping code generator 42 
receiving a chipping code clocking signal from a code clock NCO 
44. The code generator 42 may further include a prompt replica, 
not shown, for improved sampling. The chipping code may be a 
direct sequence spread spectrum code. The correlator 40 may 
operate, for example, at 50Hz. The correlator 40 provides 
inphase and quadrature, early and late, GPS correlated sample 

signals I E , Qg, I L , and Q L , communicated, for example, at IK 
Hertz into a measurement residual Kalman prefilter 46 having, 
for example, a one Hertz output. The prefilter 46 can receive 
I and Q sampled outputs from one or more correlators 40 to 

accommodate dual f^ and f-2 frequency integration using, for 
example, two respective NCOs 36 and 44 and two respective code 
generators 42 , for providing an integrated error vector to the 
single integration filter 48. The Kalman prefilter 46 computes 

the residual carrier phase errors <j> e , carrier frequency error <j> e 
nd carrier frequency rate error ty e , residual code phase error r e 

and a matrix P. The carrier phase errors and <J> e may, in 

the context of inertial navigation, respectively correspond to 

a pseudorange, pseudorange rate and pseudorange acceleration. 



-15- 



1 

2 
3 
4 
5 
6 
7 
8 
9 
10 
11 
12 
13 
14 
15 
16 
17 
18 
19 
20 
21 
22 
23 
24 
25 
26 
27 
28 



The carrier frequency error <j> e is the first time rate of change 
derivative of the carrier phase error <j) e . The carrier frequency 
rate error ^ e is the second time rate of change derivative of 
the carrier phase error The Kalman prefilter 46 also 

computes the amplitude a, code phase error r e , and a measurement 
covariance matrix P as part of a measurement vector. The 

Kalman prefilter 46 computes the measurement state covariance 
matrix P indicating the uncertainties in the measurement vector 
of errors for a number of samples m over a major cycle epoch 
time. The Kalman prefilter 46 can be used as part of a weak 
lock method using local carrier and code tracking loops for 
improved code and carrier tracking. The Kalman prefilter 46 can 
further be used in an ultratight method for GPS pseudorange 
computation tightly coupled with the correlation process to 
improve the ability to maintain tracking lock upon the received 
GPS signals . 

In the exemplar weak lock GPS navigation processor, the 

correlator outputs I E , Qj;, I L , and Q L are also sampled by a GPS 
navigation calculator 49. The Kalman prefilter 46 and GPS 
navigation calculator 49 provide system outputs 51 including 
the receiver position and velocity output P/V. The code phase 

error r e is communicated to the code generator 42 in a local 
code phase tracking loop. Concurrently, the carrier frequency 

error <j> e is fed to the carrier NCO in a carrier tracking loop. 
The navigation calculator 49 computes the best estimate of the 

navigation state vector solution defined by the position and 

velocity data P/V. The GPS navigation calculator 49 computes 
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user clock error estimates (CEE) communicated to a pseudorange 
and rate calculator 56. 

In the exemplar ultratight GPS navigation processor, IMU 
samples 52 are communicated to a 100Hz inertial navigation 
calculator 54 providing the receiver position and velocity 
outputs P/V. The Kalman prefilter measurement vector of errors 
and the state covariance matrix P are communicated to the one 
Hertz integration Kalman filter 48 having a one second major 
cycle time period between Kalman filter epochs. The integration 
Kalman filter 48 in turn repetitively provides estimates of the 
position and velocity errors in an error state vector (ESV) 
that are used to compute the position and velocity output data 
50 communicated to the control unit 14 . The IMU samples 52 are 
communicated to a 100 Hz inertial navigation calculator 54 that 
computes the best estimate of the navigation state vector 
solution defined by GPS receiver position and velocity data 
P/V. The term P/V indicates the calculated position and 
velocity data. The position and velocity data P/V is 
calculated from the IMU samples 52 and the error state vector 
ESV. At the end of each of the one second major cycles, the 
Kalman filter 48 updates the ESV and then communicates the ESV 
to the navigation calculator 54 . The IMU samples 52 are 
repetitively generated, for example, one hundred times a 
second, and include differential velocity samples Av and 
differential attitude samples AG that are communicated to the 
100Hz inertial navigation calculator 54 also receiving the 

error state vector (ESV) every second from the integration 
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Kalman filter 48. The integration Kalman filter 48 computes the 
user clock error estimates (CEE) that are communicated to the 
pseudorange and rate calculator 56. The 100 Hz inertial 
navigation calculator 54 is used for strapdown inertial 
navigation processing. 

In the exemplar weak lock navigation processor of Figure 
2A, the navigation calculator 49 provides the best estimate of 
the navigation state vector solution indicated by the position 
and velocity data (P/V) communicated to the 100 Hz delta 
pseudorange and delta pseudorange rate calculator 57 that in 
turn provides delta pseudorange rate data to the carrier NCO 36 
and provides delta pseudorange data to the code NCO 44 . The 
delta pseudorange and delta pseudorange rate calculator 57 
receives calculated GPS satellite position and velocity (P-V) 
data from a 100Hz GPS satellite position and velocity 
calculator 58 and receives the user clock error estimates 
(CEE) . The term P-V references the GPS satellite calculated 
position and velocity data. The clock error estimates (CEE) 
include a clock phase error and clock frequency error that is 
the time derivative of the clock phase error. The calculated 
delta pseudorange is equal to the change in the geometric line 
of sight range to a satellite adjusted by the clock phase error 
estimates CEE and certain signal propagation compensation 
terms, such as the nominal ionospheric and tropospheric delays 
that are not shown. The delta pseudorange is a change in the 
geometric range taken within respect to the one Hertz prefilter 
46 updates applied to the carrier NCO 36 and the code generator 
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42. The geometric range is computed from the GPS calculated 
position and velocity data P-V, and from the navigation 
calculated position and velocity data P/V. The delta 
pseudorange rate is equal to the change in the relative 
geometric velocity to the satellite from the user antenna 16 
plus the clock frequency error estimate. The GPS satellite 
position and velocity calculator 58 receives timing signals 

from a timer 60 using the reference f c 24 and receives 
demodulated ephemeris data from an ephemeris demodulator 62 

using the GPS samples 30 to compute the GPS satellite position 
and velocity data P-V. The weak lock formulation accommodates 
larger errors between the replica signal and the received 
signal by observing the errors over many data samples and 
correcting the replica signal at a much lower rate than the 
sampling rate of the correlation process inputs. 

In the exemplar ultratight navigation processor of Figure 
2B, the navigation calculator 54 provides the best estimate of 
the navigation state vector solution indicated by the position 
and velocity data (P/V) communicated to the 100 Hz pseudorange 
and pseudorange rate calculator 56 that in turn provides 
pseudorange rate data to the carrier code NCO 36 and provides 
pseudorange data to the code NCO 44. The pseudorange and 
pseudorange rate calculator 56 receives calculated GPS 
satellite position and velocity (P-V) data from a 100Hz GPS 
satellite position and velocity calculator 58 and receives the 
user clock error estimates (CEE) . The term P-V references the 
GPS satellite calculated position and velocity data. The clock 
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error estimates (CEE) include a clock phase error and clock 
frequency error that is the time derivative of the clock phase 
error. The calculated pseudorange is equal to the geometric 
line of sight range to a satellite adjusted by the clock phase 
error estimates CEE and certain signal propagation compensation 
terms, such as the nominal ionospheric and tropospheric delays 
that are not shown. The geometric range is computed from the 
GPS calculated position and velocity data P-V, and from the IMU 
calculated position and velocity data P/V. The pseudorange 
rate is equal to the relative geometric velocity to the 
satellite from the user antenna 16 plus the clock frequency 
error estimate. The GPS satellite position and velocity 
calculator 58 receives timing signals from a timer 60 using the 

reference f Q 24 and receives demodulated ephemeris data from an 
ephemeris demodulator 62 using the GPS samples 30 to compute 

the GPS satellite position and velocity data P-V. 

The GPS satellite position and velocity calculator 58 is 
used for satellite GPS position and velocity calculation. The 
pseudorange and pseudorange rate calculator 56 is used for line 
of sight pseudorange and pseudorange rate predictions . The 
delta pseudorange and delta pseudorange rate calculator 57 is 
used for indicating changes in line of sight pseudorange and 
pseudorange rate predictions. The carrier NCO 36 is used for 
generating carrier replica signals for quadrature demodulation 
of the received GPS samples 30. The code generator 42 is used 
for early and late code replica signal generation. The 
correlator 40 is used for integrate and dump correlation of the 
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I and Q carrier demodulated signals correlated with the early 
and late code replica signals. 

The correlation process is based upon the early and late 
replica code signal generation by the code generator 42 and 
upon the replica carrier generation by the carrier NCO 36. As 
the pseudorange and pseudorange rate data from the calculator 
56 are refreshed during each cycle, the replica carrier cos and 
sin signals from the carrier NCO 36 and the early and late 
replica codes signals from the code generator 42 are adjusted 
under close loop processing. The correlator 40 may operate at 
50Hz enabling conventional early and late correlation by 
providing the early and late I and Q samples to the Kalman 
prefilter 46. The Q and I samples contain amplitudinal 
informational that can be used with signal correlation 
functions to estimate the code and carrier errors . The carrier 
and code phase offsets of the replica signals as compared to 
the received GPS sample signal can be determined as offsets 
knowledge of the signal correlation function. The carrier and 
code phase offsets are used to estimate the residual errors. 

Referring to all of the Figures, and more particularly to 
Figure 3, the measurement residual Kalman prefilter 46 is based 
on a GPS signal residual estimation process that receives the 

process inputs I E , Qe, I I , and Q L that are taken as the 
correlator early (E) and late (L) , and, inphase (I) and 

quadrature (Q) samples from the correlator 40. The process 
inputs I E , Qj., I I , and Q L may be processed at a rate of lKHz. 
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The process inputs are represented collectively by the symbol 

set v=[Ie/Qe/ • Tne GPS signal residual prefilter 46 has 
four process components or algorithms, including a lKHz 

measurement residual generator 82, a lKHz state estimate 

calculator 84, a state propagator 86, and a 5Hz state error 

covariance estimator 88 . The GPS signal residual estimator 

process receives the process input y and generates the residual 

measurements § e , r e , and a. 

The measurement residual generator 82 computes a 
measurement residual r at lKHz from the process input y and the 

propagated state x. The propagated state x has five components 
where x = § e , §e, r e , a] . The propagated state x includes 

the carrier phase error state § e , the carrier frequency error 

T — 

state (|)e, the carrier frequency rate error state (j> e / code phase 
error state r e that is related to the range from the receiver to 
the satellite, and the signal amplitude state a. 

Each component of the measurement residual r depends on 

the carrier phase error state (J) e , code phase error state r e , the 
early to late replicated code offset 6 r , the broadcast data d, 

and the signal amplitude state a of the received spread 
spectrum signal 30. There is a nonlinear relationship between 

the residual components and the measurement residual r defined 

by an r matrix. 
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I e - adc(r e + 5 r ) cos <j> e 
Q e - adcfe + 5 r )sin <j> e 
I x - adc(r e - 8 r )cos$ e 
Q 1 - adc(r e - 5 r )sin <j) e 



In the residual matrix, the code correlation function c(») 
depends on the pseudorandom noise code chip width X ca as 
defined by a c(p) pseudorandom noise code equation. 



c(p) = max 



ca J 



In the c(p) pseudorandom noise code equation for a coarse 
acquisition (CA) receiver, the nominal code rate (D ca is known 
as the chipping rate where CD ca =l . 023xl0 6 chips/sec. Because the 
speed of light is C=2 . 99792458xl0 8 m/sec, the CA code chip width 
is A, ca =293.1m. The true code rate differs slightly from © ca due 
to user to satellite relative motion and atmospheric 
transmission effects. For the case where the early and late 
codes are offset from the correlation peak by 1/4 chip, 

5 r =X ca =73.26m. The broadcast message data is d=±l and has a 
nominal 20.0 millisecond period. 



A state update calculator 84 computes the state estimate x 
at a rate of lKHz . From the measurement residual r, from the 

propagated state x, and from a state error covariance estimate 
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P, the state estimate x is equal to x+K(x,P)r, where K is the 
gain. The gain K depends on components of the propagated state 

x and depends on the Riccati matrix P. The gain K is computed 
to minimize the error variance as a linear Kalman filter gain 

K(x, P)=Ph (x)V 8 t , where h is a transposed linearized 
observation matrix and V -1 is an inverse observation covariance. 
The term 8 t is a time step and may be, for example, 0.001 

second. The term h(x) is a linearized observation matrix. The 
term V is the measurement noise variance as used in an 

observation equation of the y process input. 

4 

The filter measurement yeR , where y is the process input 
including observations I E , Qe, I I , Q Lf is modeled as a nonlinear 
function of the true state x Q , which is not known, with additive 
noise v as y=G(x)+v with G(x) given by a G(x) observation matrix 
of the true state x Q . 

adc(r e + S r ) cos <|> e 
adc(r e + 8 r ) sin (|) e 
adc(r e - 8 r ) cos <J) e 
adc(r e - 8 r )sin <|> 



G(x) = 



In the G{x) observation matrix, the amplitude a, the code 
phase error r e , and the carrier phase error (|) e components can be 
a true state x 0 having a Q , r 0 , and <}> 0 components. The process 
input y is the sum of the observation matrix G(x) and the 



-24- 



1 

2 
3 
4 
5 
6 
7 
8 
9 
10 
11 
12 
13 
14 
15 
16 
17 
18 
19 
20 
21 
22 
23 
24 
25 
26 
27 
28 



additive noise v, as in an observation equation of the y process 
input . 



y = 



Ie 

Qe 

Qi 



"adc(r e + 8 r ) cos <j) e 
adc(r e + 8 r ) sin <|) e 
adc(r e - 8 r ) cos (|) e 
adc(r e - 5 r ) sin <|> e 



u 



+ 



qe 



ql 



4 

The additive disturbance veR is taken to be a random 
process with a zero-mean Gaussian distribution and the 
observation covariance V. The observation covariance V reflects 
correlation with respect to the early and late, inphase and 
quadrature correlator output signals but with the additive 

disturbance v where v=[v ie , v qe , v ±1 , v ql ] , where v is 
uncorrelated both in time and with respect to the I and Q 

correlator output signals. With 8=1/2 being the offset, 
normalized to the code wavelength, between the early and late 
correlation signals, then the observation noise covariance V is 
given by an observation expectation covariance equation. 



V = e[uu t ] = (ari) 2 



1 0 y 0 

0 1 0 y 

y 0 1 0 

0 y 0 1 



The observation expectation covariance equation includes a 
measurement noise structure matrix including the terms 0, 1 and 
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y that define a measurement noise structure. In the covariance 
expectation equation, E is a statistical expectation of the 

additive disturbance v, y=l-8, the scalar a is the signal 
amplitude is a component of the state x, and l/r\ is the signal- 
to-noise ratio. The observation covariance matrix V can be 

shown as an exemplar observation matrix V Q . 



(0 . 3f 



0.990 
0.004 
0.478 
0.218 



0.004 
1.004 
0.227 
0.488 



0.478 
0.227 
1 .001 
0.005 



0.218 
0.488 
-0.005 
1.005 



The state estimation error e is defined as e=x Q -x as in a 
state estimation error equation. 



A - Kh^ 



Kl) + G) 



In the state estimation error equation, 4 is the time 
derivative of the state error estimate. The process input 

measurement y is expressed as a function of the state 

estimation e in a y expansion equation. 



dG 

y = G(xJ + u - G(x) + — 

Ox 



e + 



+ u 
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The the expansion equation, y expanded about the 

propagated state state x. In the expansion equation, y is the 
process input from the correlator 40 and the first partial 

derivative of the observation matrix G(x) is the linearized 

observation matrix h(x) and is defined by an h(x) partial 
derivative equation. 

*W = -r- 

The h(x) partial derivative equation can be written in 
matrix form as an h(x) linearized observation matrix equation. 

- adcfe + 5 r )sin$ e 0 0 adc'fe + S r )cos<j> e dc(r; + 8 r )cos<j) e 
adcfe + S r )cos(j) e 0 0 adc'(r e + 8 r )sin(|> e dc(^ + 8 r )sin(j) e 

- adcfe - 5 r )sin<j) e 0 0 adc'fe - 8 r )cos<j) e dc(f; - 8 r )cos<j) e 
adc(r; - 5 r )cos$ e 0 0 adc'(^ - 8 r )sin<j) e dc(^ - S r )sin<j) e 



h(x) = 



In the h (x) linearized observation matrix equation, c T (r e 

+5^)and c T (r e -5y) are the code correlation sensitivities to range 
error, and can be expressed by range partial differential 



equations , 



c'(r e + 8 r ) = c(r e + 8 r ) 
c'fe - 5 r ) = A c(r e - 8 r ) 



dc(p) 



o |p| > K a 

^ o > P > -(^ ca + l) 

- rir o < p < ta a + i) 
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The state propagator 86 computes at lKHz the propagated 

state x where the propagated state components are § er $e/ <j)e/ r e 
, a. The state propagation dynamics are linear, such that, <j) e 

=4e + <Mt + 4>e8t/ <t>e=(f>e+(j)e8t i <j>e=(f>e/ r e =r e + (A, L1 /27t) S t <j)e, and a=a, 
wherein 5 t =0.001sec and 8 t is the propagation period, and 

X, L1 sO. 1903m, that is, A, L1 is the carrier wavelength. The LI band 

carrier frequency is © L1 and cd l1 =2 . 99792458x10 Hz, so that 

X L1 =C/o L1 =0. 1903m, where C is the speed of light. 

The state error covariance estimator 88 provides a state 

covariance estimate P at 5Hz preferably only using the x 

propagated state components r e and a. The state covariance 
estimator 88 may be implemented as a look-up table. Elements 

of a P scheduling look-up table, that could be a trivial look- 
up table of only one element, are computed using a state 
dynamic matrix A. A continuous time representation of the state 

dynamics is linear with x=Ax+o, with x being the time 
derivative of the state x, and with A being the state dynamics 

matrix. 





0 


1 


0 


0 


0 




0 


0 


1 


0 


0 


A = 


0 
0 
0 


0 
0 


0 
0 
0 


0 
0 
0 


0 
0 
0 
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In the state dynamics matrix A, A. L1 =0. 1903m and is the 
carrier wavelength. The additive disturbance coeR is taken to 
be a random process with a zero-mean Gaussian distribution and 

state covariance W. The state error covariance matrix P is 
provided by a continuous time Riccati equation. 

P = AP + P T A + W - Ph T V _1 hP 

In the continuous time Riccati equation, P is the time 
derivative of the Riccati state error covariance matrix P, A is 

the state dynamics matrix, V -1 in an inverse observation 

T 

covariance, h is the linearized observation matrix, h is a 
transposed linearized observation matrix, W is the state 
covariance, and P is the Riccati matrix. The measurement noise 
structure, as in the covariance expectation equation, leads to 

T -1 

a product h (x)V h(x) that does not depend on carrier phase 

error <|) e nor the carrier frequency error <j> e . The inverse 
observation covariance V 1 is given by an inverse covariance 
expectation equation. 



V- 1 = (an)" 2 



1 


0 


- Y 


0 




0 


1 


0 


- Y 


1 


- Y 


0 


1 


0 


1-Y 2 


0 


- Y 


0 


1 





In the covariance expectation equation, y-(l-5) and 5=1/2 
is the normalized early to late code offset 5 r . 
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T -1 

The product h V h is then 
equation. 

h T (x)V 1 h(x) 



defined by a scalar covariance 





0 


0 


0 


0 
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0 
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0 
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0 


0 


0 


0 


0 


0 




q ra 


0 


0 


0 


q ra 





In the scalar covariance equation, the scalar covariance 

terms q^, qyy, qy a / anc * "^aa are defined by scalar covariance 
equations . 



<J* = a2 ( C p + 2 YS C m + C m) 

q rr = a 2 (c' p 2 + 2yc' p c' m + c' m 2 ) 

q aa = ci + 2yc p c m + cl 

q ra = a L c pC' p + y(c p c' m + c ra c' p ) + c m c' m ] 



In the scalar equations Cp=c(r e +8 r ) , CjS>=c' (r e +5 r ) , c^cir^- 

8 r ) / ci =c(r e +8 r ) * The term h V h does not depend on message 
data d=±l because the message data d is a multiplication factor 
of each term in h (x) and d 2 =l . The product h T V" 1 h depends only 

T -1 

on range error r e and signal amplitude a. Because h V h does 
not depend on the harmonics sin<j) e or cos<|) e , the continuous time 
Riccati equation has steady-state solutions with respect to 

constant values for r e and a. Point design and gain trajectories 
can be determined with respect to the range error r e and the 
amplitude a* A constant error dynamics solution is determined 
by setting the range error and the signal amplitude to nominal 
values r e =0 and a=l . 
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T -1 

Because the product h V h depends only on the range error, 
that is, the code phase r e and the signal amplitude a, a steady- 
state solution of the continuous time Riccati equation can be 
expressed as an algebraic Riccati equation. 

0 = AP + P T A + W - Ph T V _1 hP 

Nominal closed-loop filter dynamics are found by setting 
the state covariance matrix W in the observation covariance 

matrix V to place the closed-loop filter poles at a(A-Kh)={- 
335,-152,-106,-8.7±4.9i) . 



W 



10 0 

0 10 3 

0 0 

0 0 

0 0 



0 
0 

10 6 

3 x 10 5 
0 



0 0 
0 0 
0 

10° 0 
0 1 



-3x10 

i5 



V= 10"% 

In the state covariance matrix W, the cross terms are 
found to improve the damping in the carrier frequency error 

rate <j> e response. The process outputs are the state estimate x 
and state error covariance estimate P with the state estimate x 

including the § e , <^ e §er r e , and a components. 



The use of the algebraic Riccati equation enables state 
error covariance estimation in real time in the Kalman 
prefilter 46. The Kalman prefilter 46 provides code phase error 
and carrier phase error signals for adjusting the code phase 
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and carrier phase for improved code phase tracking and carrier 
phase tracking. Those skilled in the art can make enhancements, 
improvements, and modifications to the invention, and these 
enhancements, improvements, and modifications may nonetheless 
fall within the spirit and scope of the following claims. 
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